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We numerically examine the nonlinear rubber elasticity of topologically constrained polymer networks. 
We propose a simple and effective model based on Graessley and Pearson's topological model (GP model) 
for describing the topological effect. The main point is to take account of a nonequilibrium effect in the 
synthesis process of the polymer network. We introduce a new parameter y to describe entropic contribu- 
tions from the entanglement of polymer loops, which may be determined from the structural characteristics 
of the sample. The model is evaluated in the light of experimental data under uniaxial and biaxial deforma- 
tions. As a result, our model exhibits uniaxial behaviors which are common to many elastomers in various 
deformation regimes such as Mooney-Rivlin's relation in small extension, stress divergence in the elon- 
gation limit and the declined stress in compression. Furthermore, it is also qualitatively consistent with 
biaxial experiments, which can be explained by few theoretical models. 

KEYWORDS: rubber elasticity, polymer network, topological effect, Mooney-Rivlin's relation, linking proba- 
bility 



1. Introduction 

The topological effect in polymer systems has been attract- 
ing considerable attention in fundamental studies as well as in 
applications. 1_5) The topologically constrained systems such 
as ring polymers have been intensively studied. 1 ' 2 ' 4-8 ) Recent 
studies 9-12) have evidenced that several familiar polyconden- 
sated polymer systems, e.g., polyesters, polycarbonates, poly 
(ether- sulf on) s and poly amides, also contain large amounts of 
cyclic polymers. Such topology of polymer chains have great 
influences on the systems. For instance, the relation between 
the length and the size of an ideal unknotted ring polymer 
shows a scaling exponent larger than that of an ideal linear 
polymer. It can be explained such that a polymer having the 
trivial knot topology behaves as though it had an excluded 
volume; that is, the constrained topology statistically results 
in expansion of the ring polymer. 13 ) It is also well known that 
the topological effect is closely related to the osmotic 14-16 ) 
and viscoelastic 17-19) properties. 

The topological constraint is expected to influence entropic 
behaviors of rubber-like materials and cause the complex- 
ity and nonlinearities of their elasticity. A rubber consists of 
polymer strands cross-linked to each other and the backbones 
form the network structure. Strands in the networks cannot 
self-intersect and intersect with others, and therefore the sys- 
tem is topologically kept unchanged under a deformation. 
Such systems demonstrate a nonlinear rubber elasticity; 20 ' 21 ) 
i.e., the uniaxial stress does not linearly depend on the strain 
for large elongation. 

The stress-strain relation in rubber-like materials is given 
by the strain energy density W, which is defined as the energy 
stored in an elastic material under a deformation: 20 ' 21 ) 
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(i) 



uniquely in terms of the invariants of Green's deformation 
tensor 



1 1 = /I 2 + + /i 2 , 
7 2 — A^A^ ~*~ $\$\ 



(2) 
(3) 



where A\ , A2 and A3 are the deformation ratios in the princi- 
pal axes, and dj are material constants. Classical theories as- 
suming an ideal chain network without molecular interactions 
gives the function W in the form 



(4) 



(5) 



Here the deformation of the material can be expressed 



W = Cio(/i - 3). 

It leads to the neo-Hookean stress 

cr = 2Ci, (/I - X' 2 ) , 

where A is the strain ratio in the stretched direction and the 
incompressibility of a sample is assumed. This stress (5) is 
correct only in the very small extension regime as shown be- 
low. 

Figure 1 shows the relation between the normalized uni- 
axial stress versus the inverse of the strain ratio (which 
we refer to as the Mooney-Rivlin plot hereinafter) obtained 
from experimental data of end-linking poly (dimethylsilox- 
ane) PDMS network. 22 ) We calculated the stress from the 
approximate functional form of W suggested by the experi- 
ment, 22 ) 

W = C 10 (h - 3) + Coi(h - 3) + C n (h - 3)(/ 2 - 3) 

+C 20 (/i - 3) 2 + C 02 (/ 2 - 3) 2 , (6) 

along with the values of Qj reported in Ref. 22) The neo- 
Hookean stress (5) obviously cannot explain the Mooney- 
Rivlin plot, giving only a flat line. In the small nonlinear ex- 
tension regime (0.5 < A~ l < 1), the stress-strain relation is 
determined almost exclusively by the first- and second-order 
terms of the function W 
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W = C w (h - 3) + Coi(/ 2 - 3). 



(7) 
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Fig. 1. The normalized uniaxial stress versus the inverse of the strain 
curves. The regimes A -1 < 1 and A~ l > 1 of the horizontal axis corre- 
spond to elongation and compression, respectively, and the case A' 1 = 1 
corresponds to the undeformed state. In the small extension regime (0.5 < 
A < 1), the curve approximately gives the straight line of Mooney-Rivlin's 
relation (8) as shown in the inset. The stress diverges in the elongation limit 
A/A c . 



It leads to Mooney-Rivlin's phenomenological equation 

(T 



21) 



- 2Cin + 2Co\A 



-l 



— -v III ' - ^ "I' 1 •> (8) 

A — A 

where Cio and Coi are positive constants. This would give 
an increasing line in the Mooney-Rivlin plot in Fig. 1. In the 
inset, we see that the stress-strain relation indeed obeys Eq. 
(8) in the small extension regime 0.5 < A~ l < 1. 

However, the second- and higher-order terms are needed 
for describing the entire deformation regime. In the large ex- 
tension regime (A~ l < 0.5), the stress <x diverges more rapidly 
with respect to A than Eq. (8) predicts and the system reaches 
the elongation limit A f A c at a smaller deformation A c than 
the ideal network does. It can be explained by adding fourth 
and fifth-order terms in Eq. (6). Further, Eq. (8) does not re- 
produce the behavior in the compression regime (A~ l > 1), 
either. The compressive stress generally gives a declining 
curve 21 ' 24) as shown in Fig. 1, contrary to Eq. (8). 

It is one of our aims in the present paper to explain 
the above nonlinearities, which correspond to higher-order 
terms of the function W, in terms of the topological effect. 
There have been several theories 20, 2 1,25_28) accounting for 
the nonlinear rubber elasticity from the point of view of the 
constrained topology of the system. However, the way that 
quenched topologies affect the networks is still not fully un- 
derstood mainly because of their complex structures. To dis- 
cuss the effect, we will extend Graessley and Pearson's topo- 
logical classification model (GP model) 27) by introducing a 
new assumption associated with the nonequilibrium process 
of the system formation. As a result, our new model qualita- 
tively reproduces the experimental facts both in extension and 
compression regimes as reported in Fig. 6 below. 

The paper is organized as follows. In section 2, we review 
two topological models: the GP model 27) and the simple net- 
work model. 28) In section 3, we explain our new model based 
on the GP model with a new assumption. In section 4, we 
numerically evaluate the uniaxial and biaxial stresses and dis- 
cuss the validity of our model. Section 5 is devoted to conclu- 
sions. 
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Fig. 2. The trivial link and some nontrivial links. The Hopf link and more 
complicated link;, i.e., topologies other than the trivial link, are called the 
nontrivial links. 



2. Topological Models 

In the present section, we review two models that partly 
explain the nonlinear rubber elasticity. The main purpose of 
the present paper is to merge the two models and to explain 
the elasticity in the entire regime. 

2. 1 Graessley and Pearson 's Topological Model 

The GP model 27) predicts the elastic contributions of poly- 
mer loops having different topologies such as entangled and 
nonentangled. The model takes ring-shaped structures, or 
loops, as structural elements of a network and employs the 
following approximations and simplifications: 

(i) Loops are randomly distributed in space with the density 
Po. 

(ii) Loop pairs are topologically classified into two cases: 
mutually entangled with the fraction P^ and not entan- 
gled with the fraction 1 - Pl- 

(iii) Loops are not entangled to themselves, i.e., they have the 
trivial knot topology. 

(iv) Distortion of loops are neglected. 

(v) Networks are affinely deformed by the external force; 
that is, the vector connecting the centers of mass of the 
two loops, r, is affinely transformed. 

(vi) The entropic contributions of loop pairs are independent 
and additive. 

In addition, the GP model assumes the existence of a charac- 
teristic distance S . In the following, we will use the radius 
of gyration of a loop as S and take the normalized vector 
£ = r/S = (x, y, z) instead of r. 

The topology of a loop pair is expressed by a link type. 
Some links are depicted in Fig. 2. Under the assumption (ii), 
the GP model distinguishes only two kinds of topologies, the 
trivial link and the nontrivial link, neglecting detailed topo- 
logical types of the nontrivial link. The topologies of the loop 
pairs are kept unchanged under a deformation; that is, the 
topological constraint condition exists. 

Here, we introduce the linking probability of a pair of 
two loops to estimate the entropy of the network. The link- 
ing probability Pl(0 is defined as the following normal- 
ized configuration number of entangled loops; Pl(0 = 

ftiink(^/(fttriv(f) + Q iink(£)X where £W£) and £1^(0 de- 
note the accessible phase space volume of a loop pair sepa- 
rated in distance £ with the trivial link and with the nontriv- 
ial link, respectively. Here we note that each loop has only 
the trivial knot topology under the assumption (iii). Thereby, 
the entropies of a loop pair having trivial and nontrivial link 
topologies are defined by s t ri v (£) = &b ln(l - Pl(0) an d 
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siink (C) = k B IiiPlC^), respectively, where k B is the Boltzmann 
constant. 

After a network reaches the equilibrium state, polymer 
loops may behave statistically as randomly distributed loops. 
Then, the assumption (iv) allows us to calculate Pl(0 both 
in deformed and undeformed states as the probability for ran- 
domly distributed undistorted loops to be entangled. It can 
be numerically obtained through the simulation of a random 
polygon (RP) or a self-avoiding polygon (SAP). 29) An SAP 
is modeled by a sequence of N line segments of unit length, 
namely Kuhn segments, and hard spherical beads of radius r^, 
which are placed at vertices to express the excluded volume. 
The beads are prevented from intersecting with one another. 
The case for r d = thus corresponds to an RP. Hence, Pl(0 
can be given by the linking probability of SAPs and RPs. 

Under the assumption (v), £ changes to = X • g if a sam- 
ple is deformed with the deformation tensor X - (A x , A y , X z ). 
The entropy change of a nontrivial link is then given as fol- 
lows: AsunkCM = *Hnk(0 " SlinkGD = *B ln(P L (f )/Pl(0). 

The entropy change of a trivial link is similarly given by 
As triv (r, = k B In [(1 - P L (D/(1 " Pl(01 The GP model 
then assumes the average entropy change per loop pair of dis- 
tance £ as follows: 



Asi oop (r , = PL(D&SKnk(D + (1 " PL(0)^ triY (0 



PL(f)ln§^ + (l-P L (^))ln 



(9) 



puo ' i-pao 

Note here that these calculations are done in the framework of 
equilibrium statistical mechanics. We will revise this equation 
in section 3 below. 

Let us consider that an external force is applied along 
the x-axis. The separation distances for undeformed and de- 
formed states are given by f 



^x 2 +y 2 +z 2 



and 



£' = \X - C\ - ^(Ax) 2 + A~ l (y 2 + z 2 ), respectively, where 
A = A i is the stretch ratio in the x-direction. (Here we 
assume the incompressibility of the system and then have 
Ai - A3 = /i~ 1/2 .) Nominal uniaxial stress is then given by 
cnoop = -(T/V)(dAs\ 00V /dX) 9 where T is the temperature. We 
thus have the following stress of a loop pair caused by the 
topological constraint: 



cr loop (£',£) = - 



k B T 
V 



Pl(0 1 - Pl(0 



pan i-plU') 
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(10) 



Under the assumption (i), one central loop is paired with 
poS 3 dxdydz pieces of other loops within a volume element 
S 3 dxdydz- Here po is the loop density defined as po = g/V, 
with the number of loops £ and the system volume V. Fur- 
thermore, under the assumption (vi), the total stress of loop 
pairs which the central loop participates can be calculated by 
summing up the contribution from individual pairs all over the 
system, such that poS 3 § y dxdydzAo-\ 00V . Then, the topologi- 
cal stress of the entire system A<xgp can be given by: 



<tgp = 



k B Tp 2 S 2 



Jo Jo Jo 



dx dy dz 
Pl(C)-Pl(0 dPiiC) 



' pana-pao) OA ' 

The pref actor of 1 /2 compensates for the over-counting of 
loop pairs. 
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Fig. 3. Linking probability Pl(£, N, rj) versus distance £ for SAPs of N = 
256 with the following seven values of the excluded volume: rj = 0.0 
(open circles); = 0.05 (closed triangles); = 0.1 (open squares); rj = 
0.15 (closed circles); = 0.2 (open triangles); rj = 0.25 (closed squares); 

= 03 (open crosses). Solid lines are fitting curves given by the formula 
(12) with the parameter values reported in Ref. 29) 



The evaluation of the GP stress requires the precise func- 
tional form of Pl(0 since Eq. (11) strongly depends on it. 27) 
The f -dependence of Pl, however, had remained obscure and 
the validity of the GP model had not been fully examined. 
In the previous study, 29 ) we numerically calculated the link- 
ing probability of off-lattice SAPs and proposed the following 
approximate formula of P^\ 



Pl&N, r d ) = exp(-< 1 ) - Cexp(-#T 2 ). 



(12) 



The fitting parameters C, a, f3, v\ and V2 depend on N and rj. 
It gave remarkably good fitting curves to the numerical results 
from the viewpoint of x 2 values, as shown in Fig. 3. 

The GP stress had not been obtained with a well fitting for- 
mula of Pl- To precisely estimate it, we here applied the well 
fitting formula (12) to Eq. (11) and obtained the result in Fig. 
4 (a). The curve reproduces the stress divergence in the elon- 
gation limit A / A c as generally observed in real networks. 
However, the small extension regime (0.5 < A~ l < 1) is 
contradictory to Mooney-Rivlin's phenomenological law (8). 
The compressive stress regime (A~ l > 1) is not consistent 
either with the experiments, 19 ' 21 ' 24) which generally show a 
decrease of the stress at A~ l increases as shown in Fig. 1. 

2.2 The Simple Network Model 

Tsurusaki 28) numerically showed that even a network with- 
out entangled loops reproduces Mooney-Rivlin's relation (8). 
We hereafter refer to it as the simple network model. The en- 
tropy change per loop pair was then given as follows: 



As ]oop (C,0 = k B ]n \ p L( /' } . 

1 - P\XO 



(13) 



The first term in the square brackets of Eq. (9), giving the 
contribution of nontrivial link pairs, was omitted here because 
the network was assumed to have no such loops. 

Further, the density of loops paired with the central loop 
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Fig. 4. The Mooney-Rivlin plot for data obtained from (a) the GP model 
and (b) the simple network model using the function P(g) of Eq. (12) with 
the fitting parameters for r j = and N = 256. 29) 



was rewritten by the following function: 

P (0 = Po-y- l e- 2 ^ 2 



(14) 



The second new term in Eq. (14) was introduced to accom- 
modate the assumption (i). It expresses that the central loop 
has few other loops around itself. In contrast, the original GP 
model sets p to be the constant po = £/V. It is noticed that 
p(f) increases with f towards po. The integral of Eq. (14) over 
the system gives £ - 1 , which is the number of loops except 
the central loop. Under these new assumptions, the uniaxial 
stress was given as follows : 28) 



k B Tp S 3 



o Jo Jo 



dx dy dz p(0 



1 



OA 
(15) 



Figure 4 (b) shows the elastic behavior of the simple net- 
work model obtained with the well fitting function of Pl(&, 
Eq. (12). As a result, the uniaxial stress demonstrated the 
Mooney-Rivlin elasticity shown in Fig. 1 for small extensions 
(0.5 < A < 1). This suggests that, surprisingly, the additional 
Coi term in Eq. (8) results from network structures without 
entangled loops; i.e., the topological interaction of trivial link 
pairs, only. In addition, the compression curve in A~ l < 1 
was qualitatively consistent with experimental facts in Fig. 1. 
However, the simple network model did not reproduce the di- 
vergence of the stress in the elongation limit A / A c . 

To summarize so far, the GP model with the well fitting 
formula of Pl(0> Eq. (12), can express the elongation limit 
A /* Ac, but not the small extension regime or the compres- 
sion regime. The simple network model with the new assump- 
tion (14), on the other hand, explained Mooney-Rivlin 's rela- 
tion in the small extension regime as well as the compression 
regime. Nevertheless, it did not reproduce the elastic behav- 



ior under large stretching for A~ l < 0.5. These results imply 
that both models cannot well represent higher-order terms of 
the phenomenological function W in Eq. (6). The main pur- 
pose of the present paper is to explain it in the entire deforma- 
tion regime by introducing a new assumption in the following 
section. By the modification, the model demonstrates good 
agreement with experiments. 

3. The GP model with a New Assumption 

We improve the GP model by introducing a new assump- 
tion in order to describe the entire deformation regime of the 
rubber elasticity. Considering possible problems of the GP 
model, the most significant would be that the fraction of en- 
tangled pairs in a system was implicitly postulated to be equal 
to the linking probability of SAPs. 

A network may be created during a highly nonequilibrium 
synthesis process. Then the network entails some loop entan- 
glements with an unknown probability, which should be dif- 
ferent from Pl(0 that reflects equilibrium properties. In or- 
der to take the nonequilibrium process into account we now 
introduce a new function /(f) as the fraction of entangled 
pairs in the system when they are separated in distance f . 
Then, the fraction of trivial link pairs is given by 1 - /(f). 
The GP model implicitly assumed that the function /(f) can 
be replaced by Pl(0 of SAPs. 27) However, such a treatment 
would be possible only if the synthesis occurred under the 
equilibrium state. The topology of loop pairs is in fact deter- 
mined through the network synthesis process, when the sys- 
tem is under a strongly nonequilibrium condition. Therefore, 
the function /(f) must, in principle, differ from Pl(0- On 
the other hand, we assume that the system is under the equi- 
librium condition during the deformation process. It corre- 
sponds to the situation where the elastic behavior of the net- 
work is observed after time enough to reach equilibrium of 
the synthesized system. Accordingly, the entropy change of a 
loop pair is given by the equilibrium function Pl(0 of SAPs, 
such as Acinic = In (Pl(<T)/^l(<D) f° r an entangled pair and 
As triv = ln(l - - P\iO)) for a not-entangled pair. 

Figure 5 illustrates the distinction between /(f) and P^. 

The contributions from trivial and nontrivial link topologies 
to the average entropy change, thereby, should be given as 
follows: 

A* loop (r , o = f(OAs link (c, o + (i - /(? ))A* riv (r , o 



i - Pl(0 



(16) 



Pl(0 v ~ " 1~Pl(0 

This is the main point of the present paper. Note again that 
each entropy change for the trivial and the nontrivial links, 
Acinic and Astdv, should still be given by the equilibrium func- 
tion Pl(0 since we consider the elasticity under the equilib- 
rium condition. 

Let us recall that the phenomenological Mooney-Rivlin law 
(8) can be explained by the simple network model without 
loop entanglements. It suggests that actual networks possi- 
bly have far fewer nontrivial linking than implied by Pl» al- 
though / and may have similar ^-dependence. Therefore, 
we introduce the following phenomenological assumption; 
the prefactor /(f) of Eq. (16) is approximately expressed by 
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Fig. 5. Linking probabilities of a loop pair (a) in the nonequilibrium syn- 
thesis process, /(£) and (b) in the equilibrium state, Pl(0- The function 
/(£) indicates the probability that two loops of distance £ are entangled to 
each other through the nonequilibrium synthesis process of a network. On 
the other hand, P\££) is defined as the linking probability of two equilib- 
rium SAPs. After the network reaches the equilibrium state, the number of 
the configurations of the loop pair would be given by Pl(0- Therefore, we 
express the entropies of a loop pair in the undeformed and deformed sys- 
tems as ^lnPL(^) and mPi^X respectively. We thereby obtain the 
mean entropy change per loop pair by deformation as in Eq. (16). 



the form 



f(0 = 7 P L (0- 



(17) 



Hence, the parameter y is a quantity giving the entropic con- 
tribution of the entangled loops. It should depend on the struc- 
tural characteristics of a network. The above consideration 
leads to the following average entropy change per loop pair: 

As loov (C,0 



r P L (f)ln5^ + (l- r P L (f))ln 



Pl(0 

The stress contributed by a loop pair is 



1 - Pl(0 



(18) 



dA ' 



(19) 



^ioop(r^)- v [y- PL(n i-p L(r) 

We thus have the uniaxial stress of our new model, with the 
position-dependent density p(£) in Eq. (14), as follows; 

k B Tp S 3 



0"new — 



/^»00 /^»00 /^»00 

I I I dxdydzpiO 
Jo Jo Jo 



P L (£')(1 - P L ((>)) dA 



(20) 



We will also use the well fitting formula (12) for the function 



Fig. 6. The Mooney-Rivlin plot of the result of our new model (20) with 
7 = 0.1. Here we used Eq. (12) in the stress calculation with the fitting 
parameters for = 0.0 and N = 256. 



4. Results and Discussion 

4. 1 Uniaxial Elasticities in Extension and Compression 
Figure 6 shows the uniaxial stress of the new model <x new 

with y = 0.1 in the Mooney-Rivlin plot. The result is in good 
agreement with the typical experimental data of the nonlinear 
rubber elasticity throughout the entire deformation regime. 
Our new model qualitatively reproduces elastic behavior with 
higher-order terms of Eq. (6) under extension and compres- 
sion. 

The best value of the fitting parameter y can vary depend- 
ing on the length N and the excluded volume of a loop. In 
the case N - 256 and = 0.0, the stress successfully fits the 
empirical form as shown in Fig. 6 within the values of roughly 
0.01 < y < 0.1. The result suggests that the nontrivial link- 
ing rarely occurs in the synthesis process of a real network, 
contrary to the expectation implied by P^. As y increases, the 
Mooney-Rivlin curve approaches that of the GP model shown 
in Fig. 4 (a), since the number of the nontrivial links increases. 

The choice of the parameter y seems to be a feasible task 
if the network structure can be determined. It may be, how- 
ever, exceedingly difficult for real networks because their mi- 
croscopic topological structures are complex and not directly 
analyzable. 

4.2 Elastic Contributions from Trivial and Non-Trivial 
Links 

Both the simple network model and our new model, with 
a small value of y, show the stress-strain relation follow- 
ing Mooney-Rivlin 's theory in the small extension regime. 
In contrast, the stress divergence under large deformation 
can be explained only from systems containing the entangle- 
ment; such as, the original GP model and our new model with 
7^0. These results imply that the trivial and nontrivial link 
topologies govern the system elasticity in different deforma- 
tion regimes. Let us now focus on this topic. 

Figure 7 shows the contribution of trivial and nontrivial 
links to the stress. Their elastic behaviors are considerably dif- 
ferent. Under small deformations (0.5 < AT 1 < 1), the change 
of the normalized stress which is caused by the trivial link ex- 
ceeds that due to the nontrivial link. Therefore, the trivial link 
topology dominates the system elasticity in the small defor- 
mation regime and its stress causes Mooney-Rivlin 's relation. 
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Fig. 7. The Mooney-Rivlin plots of the stresses of nontrivial and trivial link 
topologies. Here we set = 0.0, N = 256 and y = 0.1. 
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Fig. 8. Topological interaction between loops having different topologies: 
(a) the trivial link and (b) the nontrivial link. 



Let us explain this behavior in terms of Eq. (19). In Eq. 
(19), the elastic contributions from nontrivial and trivial links 
are determined by P\XC)- Under small extensions, loops may 
be closely distributed to each other and their linking probabil- 
ity Pl(C) should have a large average as seen in the numerical 
results in Fig. 3 for small £. In this case, the stress of the trivial 
links drastically changes with X since the /i-derivative of the 
second term in Eq. (19) is proportional to (1 - P(£'))~ 2 and 
takes a large value when Pl(C) ~ U whereas, the stress of 
the nontrivial links slightly changes because the /1-derivative 
of the first term in Eq. (19) is proportional to P(£')~ 2 . 

In the large extension regime (A~ l < 0.5), on the other 
hand, the contribution of the trivial links vanishes as seen in 
Fig. 7, and thus the system stress is affected principally by the 
nontrivial links. The stress of the nontrivial links diverges in 
the deformation limit A c as shown in Fig. 7. This behavior can 
be also explained in terms of Eq. (19). The linking probabil- 
ity Pl(0 for small continuously decreases as £ increases 
as shown in Fig. 3. It makes the first term of Eq. (19) domi- 
nant under large deformations since the loops tend to be very 
separated from each other in such a situation. 

Figure 8 (a) schematically shows how the topological in- 
teraction between the loops of a trivial link pair contributes to 
the stress. If a trivial link pair is close to each other in the di- 
rection perpendicular to a deformation, the entropic restoring 
force acts on the trivial link, since the loops become closer 




0.8 1.0 



Fig. 9. The Mooney-Rivlin plots of our new model with the excluded vol- 
ume rd = 0.0 (top), 0.05 (middle) and 0.1 (bottom) for the Kuhn segment 
of an unit length. We set N = 256 and y = 0. 1 . 



and interrupt each other's random configurations under the 
deformation. This restoring force contributes to the stress in 
the perpendicular direction because the system is assumed 
to be imcompressible. The interaction between such a trivial 
link pair should be the origin of the Mooney-Rivlin stress ex- 
pressed as Eq. (8). Figure 8 (b), whereas, demonstrates that an 
entangled loop pair is strongly attracted to each other under a 
deformation. Hence, the system stress diverges in the large de- 
formation regime mainly because of the nontrivial link topo- 
logical interaction. 

4. 3 Swollen Effect of the New Model 

The above considerations are for small r^. In contrast, 
highly swollen networks have a large excluded volume and 
they are strongly influenced by the existence of solvents. In- 
deed, it is frequently reported that the additional term Qi on 
Mooney-Rivlin 's relation (8) decreases with the volume frac- 
tion of polymers. 2 ^ As demonstrated in Fig. 9, our new model 
is consistent with such experimental facts; the gradient of the 
Mooney-Rivlin curve decreases with r<i increases in the small 
deformation regime (0.5 < X < 1). We can understand the 
swelling effect from the Yd -dependence of Pl- The excluded 
volume Yd drastically affects P\XO m Fig- 3. When > 0.15, 
PhiO takes rather small values even for f = 0. Hence, the 
topological effect of the trivial link topology should practi- 
cally vanish since the second term of Eq. (19) is relatively 
negligible over all the deformation regimes. 

4.4 Biaxial Stress-Strain Relation 

Let us finally consider the sample in biaxial extension il- 
lustrated in Fig. 10. The results for biaxial deformation are 
plotted in Fig. 11 (a) and (b). Our new model exhibits the 
stress-stain curves that are qualitatively consistent with an ex- 
perimental result. 22 ^ According to the Ref., 23) there are few 
theoretical models that reproduce the biaxial experiment. The 
present result for biaxial deformation is a remarkable advan- 
tage of our new model. 

5. Conclusions 

We have improved the GP model as follows. We have pro- 
posed that the probability of entangled loops, /(f), must be, in 
principle, different from the linking probability of equilibrium 
SAPs, Pl(0> since the loop topologies are determined by the 
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Fig. 10. The sample is elongated in the v-direction with the nominal strain 
ratio A2 when the strain ratio in the x-direction is fixed at X\. Nominal 
stresses along the x- and ^-directions are denoted by cr\ and <X2, respec- 
tively. 
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Fig. 11. Nominal stresses along the x- and v-directions, (a) <t\ and (b) o~2, 
respectively, calculated from the new model under combinations of the 
deformation ratios A\ and A2. Here the adjustable parameter 7 is set to 
0.1 and the parameter values C, a, ft, v\ and V2 in the case of rj = 0.0 
and N = 256 are used. Each solid line stands for the A2 -dependence of 
the stress under a fixed deformation A\ . The bottom line corresponds to 
A\ = 1 . 1 , the top line corresponds to A\ = 1.9 and the intervals are 0. 1 . This 
graph shows the sum of the topological stress of our new model and the 
linear Hookean stress. We determined the coefficient C\$ of the Hookean 
term in Eq. (6) so as to agree with the biaxial experimental data. 22) 



nonequilibrium process of network forming, when the poly- 
mers are not distributed as SAPs. Moreover, we have intro- 
duced a new phenomenological parameter y in order to eval- 
uate f(0 such that f{0 = yP L (0- 

The model qualitatively well represents nonlinear rub- 
ber elasticity both in extension and compression regime for 
0.01 < y < 0.1. Here we emphasize that our model also ex- 



plains the biaxial experiments, while few models in the pre- 
vious studies can reproduce them. We thereby claim that the 
present model gives an essential step to quantitative analysis 
of the nonlinear rubber elasticity. 

The primary advantage of the present model is that it has 
few adjustable parameters. If the parameter y is given, the 
model stress is calculated by the structural and topological 
features of a network without any other adjustable parame- 
ters. Further study is required to clarify the appropriate way 
of determining the parameter y. 
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